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1.0  SUMMARY 


The  spatial  stability  of  fourteen  Ealkner-Skan  similarity  profiles  for 

the  range  P  —  0.1988  (separation)  to  0  =  1.0  (stagnation  flow)  are 

considered.  The  Orr-Sommerfeld  equation  that  governs  the  amplitude  of 

Tollmien-Schlichting  disturbances  superimposed  on  flows  having  such  profiles 

is  solved  by  a  step-by-step  integration  method.  The  rapid  error  buildup 

usually  inherent  in  such  methods  is  eliminated  by  the  use  of  the  Gram-Schmidt 

orthogonalization  procedure.  Eigenvalue  solutions  for  values  of  Reynolds 

<3 

number  Rg*  up  to  10^  were  successfully  obtained  for  all  the  cases  considered. 
Curves  of  constant  spatial  amplification  are  presented  on  the  dimensionless 
frequency  vs.  Reynolds  number  diagram  for  the  various  0*  s.  The  corresponding 
temporal  amplification  rates  are  presented  in  the  form  of  curves  of  constant 
c^  on  the  dimensionless  wave  number  vs.  Reynolds  number  diagram  (curves  of 
constant  phase  speed  cy  are  also  included).  The  critical  Reynolds  number, 
which  is  presented  as  a  function  of  0,  is  found  to  be  in  fair  agreement 
with  Pretsch's  results.  In  the  present  study,  however,  the  separation  profile 
is  found  to  possess  a  finite  critical  Reynolds  number.  Eigenfunction  solutions 
obtained  from  the  present  analysis  are  in  excellent  agreement  with  the  calcula¬ 
tions  of  Radbill  and  Kaplan  and  also  with  the  experimental  data  of  Schubauer 
and  Skramstad. 
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4.0  NOMENCLATURE 
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coefficient'  of  planar  fit  used  to  locate  zeroes  of  secular  equation 
(see  equation  (l8)) 

Gram- Schmidt  orthogonalization  constant  defined  by  equation  (ll) 
dimensionless  complex  velocity  of  Tollmien-Schlichting  disturbance, 

cA 

dimensionless  group  velocity,  c  /U 

6*  e 

dimensionless  stream  function  of  mean  flow 
characteristic  function  defined  by  equation  (6) 

Gram-Schmidt  normalization  factor  for  at  y^ 

Reynolds  number  based  on  boundary- layer  thickness,  U^fc/v 

Reynolds  number  based  on  displacement  thickness,  U  6*/v 

denotes  spatial  values 

dimensionless  time  parameter,  t^U^/6 

denotes  temporal  values 

root -mean -square  value  of  longitudinal  component  of  disturbance 
velocity 

dimensionless  longitudinal  velocity  component  of  the  mean  flow, 

U*/Ue 

longitudinal  velocity  component  in  the  external  free  stream 

dimensionless  distance  along  body  surface,  x^/6 

dimensionless  normal  to  surface  coordinate,  y.^/6 

normal  admittance  defined  by  equation  (15) 

dimensional  wave  number  associated  with  Tollmien-Schlichting 
disturbance 

dimensional  wave  number,  06 

Hartree  0,  dimensionless  pressure -parameter 

latent  root,  [cP  +  id®(l  —  c)]"^,  of  the  reduced  Orr-Sommerfeld 
equation  in  outer  region 
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displacement  thickness 

kinematic -viscosity  coefficient 

2  2  2  2 
Laplace  operator,  d  /dx  +  d  /dy 

dimensionless  coordinate  normal  to  surface 

disturbance  amplitude  function 

vector  function  space  defined  by  equation  (9) 

eigenfunction  solution 

disturbance  stream  function 

stream  function  of  the  mean  flow 

total  stream  function 


vorticity 

dimensionless  frequency, 


constraint  function  defined  by  equation  (30  ) 


SUPERSCRIPTS 

differentiation  with  respect  to  coordinate  normal  to  surface 
denotes  normalized  values 


SUBSCRIPTS 

imaginary  part  of  complex  quantity 
real  part  of  complex  quantity 
value  of  k-th  integration  step 
wall  value 

where  there  may  be  doubt,  dimensional  quantities  are  indicated 

by  (  )* 
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5.0  INTRODUCTION 


The  stability  of  a  family  of  incompressible  boundary-layer  profiles  with 
respect  to  disturbances  of  the  Tollmien-Schlichting  type  has  been  extensively 
treated  by  only  two  people,  Schlichting^  and  Pretsch^.  Schlichting 
analyzed  the  Bohlhausen  family  and  Pretsch  the  Falkner-Skan  family.  The  pro¬ 
files  of  the  Pohlhausen  family  are  convenient  for  analysis  but,  of  course, 
are  only  approximations  of  the  exact  solutions.  The  profiles  of  the  Falkner- 
Skan  family,  which  are  exact  solutions,  are  considered  in  this  report.  The 
velocity  profiles  that  occur  in  this  family,  such  as  0  =  0  (flat-plate  flow), 

0  =  —  0.1Q88  (separation),  and  0  =  1.0  (stagnation  flow),  are  of  practical 
as  well  as  theoretical  interest.  It  was  for  a  flow  over  a  flat  plate  that 
Schubauer  and  Skramstad  first  demonstrated  the  occurrence  of  the  theoretically 
predicted  regions  of  instability  as  the  first  step  in  the  transition  process 
from  laminar  to  turbulent  flow.  The  connection  between  instability  and 
transition  is  still  largely  unresolved.  Nonetheless,  a  semiempirical  method 
for  predicting  transition  developed  independently  by  Staich^  and  Van  Ingen ^ 
enjoyed  remarkable  success.  The  theoretical  analyses  for  all  these  important 
papers  were  based  on  the  stability  characteristics  of  the  Falkner-Skan  profiles. 

The  stability  of  Falkner-Skan  profiles  is  reconsidered  here  because  the 
earlier  works  by  both  Schlichting  and  Pretsch  involve  approximate  representa¬ 
tion  of  the  mean  flow  and  rely  on  asymptotic  methods  to  obtain  solutions  of 
the  stability  equation.  It  has  since  been  found  that  asymptotic  methods, 

although  useful  in  establishing  conditions  for  stability  at  high  Reynolds 

(5) 

numbers,  often  lack  accuracy  at  lower  Reynolds  numbers.  Moreover,  because 
of  the  presence  of  the  second  derivative  of  the  mean  velocity  in  the  stability 
equation,  it  is  important  to  use  accurate  velocity  profiles  in  performing 
stability  calculations.  Hie  different  results  obtained  by  Schlichting  and 
Pretsch  for  a  flow  in  the  presence  of  an  adverse  pressure  gradient  serve  to 
demonstrate  this  point.  In  this  situation,  the  upper  branch  of  the  neutral 
curve  approaches  an  asymptote,  ar6*  =  constant,  as  R^*  —  ».  For  the 
separating  profile  (0  =  —  0.1988)  the  difference  in  the  numerical  values 
obtained  was  considerable.  Schlichting  obtained  0.1*50,  whereas  Pretsch 
obtained  1.07.  This  difference  is  primarily  due  to  the  assumed  form  of  the 
mean  velocity  profile.  Schlichting  based  his  solution  on  a  polynomial 
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approximation  of  the  velocity  profile.  Pretsch  based  his  on  the  more  accurate 
profile  calculated  by  Hartree^.  With  the  exception  of  a  few  cases  (mainly 
^  =  0  for  a  limited  range  of  Reynolds  number),  the  mathematical  difficulties 
encountered  in  the  stability  analysis  have  discouraged  further  attempts  to 
improve  the  stability  charts  produced  by  these  early  investigators.  The 
usefulness  of  the  charts  and  a  desire  to  devise  the  mathematical  techniques 
necessary  to  treat  this  and  similar  problems  efficiently  prompted  a  further 
investigation  into  the  stability  properties  of  these  profiles. 

Because  of  the  limitations  of  asymptotic  analytic  methods,  there  has  been 
a  considerable  effort  devoted  in  recent  years  to  the  development  of  numerical 
methods.  These  inherently  afford  good  accuracy  at  low  Reynolds  numbers  but 
encounter  numerical  difficulties  at  high  values.  A  notable  advance  was  made 
by  Landahl  and  Kaplan^  with  their  ’’purification  scheme”.  It  cleverly 
circumvents  the  rapid  error  buildup  that  usually  occurs  in  the  step-by-step 
integration  of  the  Orr-Sommerfeld  equation.  This  equation  is  commonly  accepted 
as  the  basis  for  stability  investigations  in  two-dimensional  incompressible 
isothermal  flows.  The  purification  scheme  was  found  adequate  to  deal  with 
cases  in  which  the  Reynolds  number  based  on  boundary-layer  thickness  did  not 
exceed  10^.  In  many  practical  cases,  critical  Reynolds  numbers  far  greater 
than  10^  must  be  considered.  For  example,  in  flows  with  adverse  pressure 
gradients,  the  asymptotic  behavior  of  the  neutral  curve  is  established  at 
very  large  Reynolds  numbers.  The  Reynolds  number  restriction  thus  imposes 
a  severe  limitation  on  the  class  of  problems  that  can  be  investigated  by 
Landahl  and  Kaplan’s  method.  Nonetheless,  their  basic  approach  is  sound  and 
the  more  precise  mathematical  formulation  of  the  ” purification”  procedure 
adopted  in  the  present  analysis  has  proved  successful  in  eliminating  the 
shortcoming. 

In  the  present  investigation  the  purification  scheme  of  Landahl  and  Kaplan 
is  replaced  by  the  Gram-Schmidt  orthogonal ization  procedure This  innova¬ 
tion  has  been  successfully  applied  to  the  solution  of  the  Orr-Sommerfeld 
equation  for  values  of  Reynolds  numbers  R^*  up  to  10  .  The  velocity  profiles 
used  to  perform  the  stability  calculations  arc*  obtained  from  a  direct  numerical 
integration  of  the  boundary -layer  equation. 
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6.0  ANALYSIS 


6.1  The  Orr - Somme rf eld  Equation 


The  dynamic  equation  for  an  incompressible,  non-heat-conducting,  viscous 
fluid  in  two-dimensional  parallel  flow  written  in  terms  of  the  vorticity 
£(t,x,y)  and  stream  function  ty(t,x,y)  is 


where 


aj  cl  _  c*£ 
dx  dy  dy 


& 

dx 


(la) 


t,  =  A* 


(lb) 


The  above  equations  have  been  made  dimensionless  by  normalizing  all  lengths 
by  the  reference  length  6,  all  velocities  with  respect  to  U^,  and  all 
times  by  b/U^.  F°r  the  purpose  of  linearization,  a  small  disturbance 
i(x,y,t)  is  superimposed  on  the  mean  flow  ;  i.e., 


*  =  40(y)  +  +(x,y,t) 


(2) 


Introduction  of  the  Tollmien-Schlichting  wave 

ia  (x-ct) 

t  =  <p(y)e  (3) 


into  equation  (2)  and  substituting  the  resulting  stream  function  into  equation 
(la)  leads,  after  linearization,  to  the  familiar  Orr-Sommerfeld  equation  (here¬ 
after  referred  to  as  the  0-S  equation), 

(u  -  c)(q>"  -  C^<p)  -  U"<P  =  -  c^R  (*P""  ~  20^’f'+  ai<p)  (**) 

Henceforth  subscript  1  on  0t  is  deleted  for  convenience.  The  boundary 
conditions  assumed  for  the  present  problem  are  homogeneous.  Thus  the  solution 
of  the  0-S  equation  constitutes  an  eigenvalue  problem  leading  to  the  solution 
of  a  secular  equation  of  the  form  G(R,  a,  c)  =0. 
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6.2  Spatial  Amplification  of  Disturbance 

Traditionally,  in  the  formulation  of  the  eigenvalue  problem  disturbances 
are  considered  as  amplified  or  damped  in  time;  a  is  taken  to  be  real 
(a  =  ar)  and  the  frequency  o>  =  Oc  is  taken  to  be  complex  (a)  =  a)r  +  ia\). 
The  temporal  point  of  view  is  of  course  ideal  for  such  analysis.  This  becomes 
obvious  when  we  note  that  in  the  governing  equations  of  motion,  the  temporal 
derivatives  are  of  first  order  and  have  a  coefficient  of  unity,  whe  *eas  the 
other  derivatives,  all  of  the  spatial  variety,  are  of  higher  order  and  appear 
in  more  complicated  forms.  The  temporal  viewpoint  has  produced  theoretical 
predictions  that  are,  in  some  cases,  in  fair  agreement  with  experimental 
observations.  Much  of  the  discussion  that  follows  in  this  section  may  be 
found  in  the  current  literature  (e.g.,  reference  7),  but  will  be  presented 
here  for  the  sake  of  completeness. 


Physical  observations  of  the  breakdown  of  laminar  flows  in  a  boundary 

layer  have  revealed  the  presence  of  combined  spatial  and  temporal  growth  of 

(9) 

small  wavelike  disturbances  preceding  transition,  e.g.,  Brown  .  On  the 
other  hand,  the  forced-oscillation  experiments  of  Schubauer  and  Skramstad^10^ 
and  KLebanoff,  Tidstrom,  and  Sargent give  illustrations  of  a  purely  spatial 
growth  of  disturbances.  These  observations  necessitate  a  closer  look  into 
the  adequacy  of  the  description  of  the  growth  of  small  disturbances  as  a 
process  of  temporal  amplification.  In  this  paper  the  eigenvalue  problem  for 
the  case  of  spatial  amplification,  ^^complex*  ^real^  =  considered# 

Two  different  approaches  for  bridging  the  difference  between  the  theoretical 
and  the  observed  descriptions  of  the  phenomenon  will  be  discussed.  The  first 


is  descriued  in  a  note  by  Gaster' 


that  derives  a  relation  between  the 


spatial  and  temporal  amplification  rates  under  the  assumption  that  these  rates 

2 

are  small.  The  transformation  states  that  to  the  order  of  o>.  ,  where  0). 

irrr  lm 

is  the  maximum  value  of  ox  at  a  given  Reynolds  number,  one  may  write 


ac.(T) 


a.  (s)  - - 

1  '  C 


+  0(cp 
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(5b) 


a.(s)  =  a(T) 

<o.(s)  =  “.(t) 

where  the  disturbance  group  velocity  c  is  given  by 

o 


c 

g 


c  +  a 
r 


*cr 

§5" 


5(acr) 

da 


(5c) 


(5d) 


and  S  and  T  denote  spatial  and  temporal  systems,  respectively.  Note  that 
the  group  velocity  of  a  dispersive  wave  system  is  the  speed  at  which  the  energy 
of  the  disturbance  is  convected.  For  fixed -frequency  disturbances,  one  must 
recognize  that  although  acf  is  fixed,  its  variation  with  wave  number  does 
not  vanish  and  that  energy  may  still  be  carried  by  a  disturbance  at  a 
velocity  different  from  its  phase  speed.  The  method  is  based  on  the 
linearization  of  the  relationship  between  a  and  c  in  any  given  neighbor¬ 
hood.  This  method  of  relating  spatial  and  temporal  amplification  rates  was 
used  by  Landahl  and  Kaplan^)  in  studying  the  stability  of  laminar  flow  over 
compliant  boundaries.  Recently,  however,  Betchov  and  Criminale^1^  showed 
that  the  linear  relationship  between  Ot  and  c  assumed  by  Gaster  is  not 
always  valid  in  regions  away  from  the  neutral  curve.  In  fact,  there  may  be 
regions  for  which  the  temporal  criterion  does  not  predict  instability  although 
the  spatial  criterion  does. 

The  other  approach  to  spatial  amplification  is  of  course  the  direct  attack 
on  the  characteristic  equation, 


G(a,u>,R)  =  0  (6) 

by  regarding  a  as  a  complex  parameter  and  considering  only  real  values  for 
thus  restricting  c  to  be  proportional  to  the  complex  conjugate  of  a. 

In  this  approach,  which  is  used  in  the  present  investigation,  there  is  no 
restriction  on  the  magnitudes  of  a  and  c,  as  there  is  with  Gaster' s 
transformation. 
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It  must  be  noted  that  a  complex  a  introduces  amplifications  of  the  form 

■Oc  x 

e  i  .  A  term  of  this  form  causes  exponential  growth  of  the  disturbances  for 
either  positive  or  negative  x  for  the  case  of  nonzero  0^.  Thus,  to  develop 
a  stability  criterion,  one  must  determine  whether  the  disturbances  grow  as 
they  travel  in  the  upstream  or  downstream  direction.  Establishing  a  criterion 
for  spatial  stability  requires  knowledge  of  the  direction  in  which  the  wave 
is  amplified  and  of  the  direction  in  which  energy  is  carried  by  the  disturbance, 
as  determined  by  the  sign  of  the  group  velocity.  The  spatial  criterion  for 
stability  is  presented  in  table  I.  The  temporal  criterion  for  stability  is 
also  presented  for  comparison.  The  use  of  the  group  velocity  in  the  stability 


Table  I 

SPATIAL  AND  TEMPORAL  CRITERION  OF  STABILITY  ACCORDING  TO 
LAMDAHL  AMD  KAPLAN^ 


Temporal  Growth 

Spatial  Growth 

Stability 

c.  >  0 

a.  <  0  c  >0 

Unstable 

l 

i  g 

a.  >  0  c  <0 
i  g 

c.  =  0 

l 

o 

ii 

«H 

Neutral 

c.  <  0 

a.  >  o  c  >o 

Stable 

l 

i  g 

a.  <  o  c  <o 
i  g 

criterion  for  spatial  growth  of  disturbances  is  an  effective  application  of 
a  radiation  condition.  Hence,  disturbances  are  unstable  as  they  convect 
energy  in  the  same  direction  as  that  of  their  spatial  growth.  By  means  of 
Caster’s  transformations,  equation  (5),  it  is  seen  from  table  I  that  the  two 
criteria,  the  spatial  and  the  temporal,  are  consistent.  The  present  paper 
reports  on  the  spatial  stability  of  Falkner-Skan  similarity  profiles  for 
which  the  characteristic  function  is  solved  in  terms  of  real  values  of  o> 
and  complex  values  of  CL, 
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6.3  Solution  of  the  Orr- Somme rf eld  Equation 


When  the  coefficients  of  the  0-S  equation  vary  in  the  range  of  interest, 
as  in  the  present  case,  no  closed-form  solution  exists.  In  the  region  outside 
the  boundary  layer,  the  mean  flow  is  constant  and  the  0-S  equation  reduces  to 
the  following  equation  with  constant  coefficients: 

<jP-  (c?  +  y2)<p  +  =  o  (7) 

where  T2  =  oP  +  io®(l  -  c).  This  differential  equation  has  four  solutions 
expressible  in  terms  of  exponential  functions,  as  follows 


e 


-ay 


<p2  = 


e 


ay 


*5 


e 


n 


(8) 


where  the  real  parts  of  a  and  y  are  greater  than  zero.  If  the  solutions 
are  to  remain  bounded  as  y  -» «>,  then  cp^  and  <p^  are  the  only  admissible 
solutions.  Since  cp^  and  cp^  satisfy  the  0-S  equation  at  the  outer  edge  of 
the  boundary  layer  (i.e.,  y  =  l),  they  are  used  to  specify  initial  conditions 
for  the  corresponding  solutions  within  the  boundary  layer  itself.  A  Runge- 
Kutta  method  for  a  fourth-order  differential  equation  is  used  to  integrate  the 
0-S  equation  with  the  parameters  R  and  Oc  fixed.  A  linear  combination 
of  these  two  solutions  is  then  formed  to  satisfy  a  boundary  constraint  at  the 
wall.  A  numerical  search  procedure  is  then  used  whereby  the  value  of  the 
parameter  a  is  found  for  which  the  boundary  conditions  at  the  wall  are 
satisfied.  The  particular  values  of  R,  a,  and  c  for  which  this  is  true 
are  referred  to  as  eigenvalues.  The  problem  is  thus  solved,  at  least  in 
principle.  In  practice,  however,  two  solutions  that  are  linearly  independent 
at  the  beginning  of  the  integration  are  observed  to  become  linearly  dependent 
within  the  region  of  interest.  This  difficulty  can  be  removed  by  using 
multiple -precis ion  arithmetic,  which  is  rather  costly  in  terms  of  programmer 
and  computer  time. 
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6.4  Gram-Schmidt  Orthogonalization  Procedure 


Initially,  at  least,  the  two  solutions  cp1  and  q>^  are  linearly 
independent.  To  insure  that  they  remain  so,  the  following  function  space  is 
defined : 

♦(y)  =  (<p,  tf',  <p"')  (9) 

A  given  vector  in  this  complex  function  space  completely  defines  a  solution 
of  the  governing  equation.  Two  solutions  or,  equivalently,  two  vectors  in 
this  function  space  that  are  linearly  independent  (noncollinear)  may  be 
linearly  combined  to  construct  another  solution  (vector)  that  is  entirely  free 
of  any  one  of  the  generating  solutions  (vector).  Because  of  round-off,  4^ 
has  a  small  parasitic  error  proportional  to  $y  To  remove  the  parasitic  error 
from  the  integrated  solution,  an  auxiliary  function  ^  is  constructed  such 
that  no  component  of  is  contained  in  it.  This  is  accomplished  by  normal¬ 
izing  the  integrated  vaJues  of  at  the  end  of  each  step  of  the  integration 

ill  constructing,  at  some  y^,  the  function  as  follows: 


h(V  -  W  -  aW 


where  ^  is  the  normalized  value  and  A  is  an  as  yet  undetermined  constant. 
Since  it  is  desired  to  remove  any  presence  of  the  solution  from  the 
function  is  constructed  to  be  orthogonal  to  4>y  Hence, 


a  =  ^(yk)  •  4>5(yk) 


=  ♦jfy)  -  j,t'1(yk)  *  *5(yk)|  <t>5(y) 


Since  ^  and  ^  represent  solutions  of  the  0-S  equation,  the  purified 
solution  Ij^  is  also  a  solution.  The  procedure  may  be  extended  to  generate 
additional  independent  solutions  as  desired.  The  algorithm  is  known  as  the 
Gram-Schmidt  orthogonalization  procedure. 
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The  purification  of  the  integrated  solution  for  ^  may  be  applied  at 
any  arbitrary  step  of  the  integration,  provided  the  solution  has  not  become 
completely  dominated  by  the  parasitic  error.  .In  the  present  study  the  ortho- 
gonalization  procedure  is  applied  at  each  step  of  the  integration.  No  inves¬ 
tigation  was  made  to  determine  the  effect  of  the  frequency  of  orthonormalization 
on  the  results.  The  above  method  is  used  in  conjunction  with  the  integration 
(IBM  709*0  of  the  0-S  equation  for  a  number  of  Falkner-Skan  velocity  profiles. 
Single -precision  arithmetic  was  adequate  for  all  the  necessary  computations  in 
the  stability  analysis. 

6.5  Numerical  Determination  of  the  Eigenvalue  Criterion 

To  solve  the  eigenvalue  problem  numerically,  one  does  not  need  to  know 
the  functions  $1(y)  and  $^(y),  but  only  their  values  at  the  wall.  A 
solution  is  constructed  as  follows: 

«  -  (13) 

where  B  is  selected  from  the  requirement  that 

c*^  +  W>w  =  0  (14) 

and  D  is  some  arbitrary  constant.  The  boundary-layer  admittance  Z,  which 
is  proportional  to  <j>w,  is  substituted  for  the  boundary  condition  =  0 


-C^R* 

Z(R,  oc,  a)  - - (15) 

•  -  a2*' 

w  w 


The  eigenvalue  criterion  is  then  the  condition  that  Z  vanishes.  When  this 
occurs, 


\  =  K  =  0  (l6) 

and  the  boundary  conditions  at  the  wall  are  satisfied.  The  quantity  Z  is,  in 
general,  complex.  Hence,  equation  (15)  may  be  expressed  separately  in  terms  of 
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its  real  and  imaginary  components,  thus  yielding  two  equations  in  two  unknowns. 
There  are  only  two  unknowns,  since  R  and  ac  are  held  fixed  and  a  is  taken 
to  be  a  complex  number.  Thus  equation  (15)  may  be  written  in  the  functional 
form 


zr  =  zr(ar,  Ojj 

zi  =  zi(“r>  ai)  (17) 

The  application  of  a  standard  iteration  method  to  solve  these  equations 
is  convenient.  However,  because  of  the  acute  sensitivity  of  the  topology  of 
Z  to  small  changes  in  a,  such  a  scheme  is  inadequate  to  meet  the  general 
requirements  of  the  problem.  The  use  of  a  LaGrange-type  interpolation  such  as 
that  employed  by  Landahl  and  Kaplan  proved  unsatisfactory  for  the  Reynolds 
number  range  considered  in  this  study.  The  procedure  now  used  is  to  fit  a 
plane  to  three  successive  trial  points  of  the  form 

a  =  aQ  +  &jZr  +  a2Zi  (18) 

where  the  a*s  are,  in  general,  complex.  The  first  two  trials  are  based  on 
either  estimates  or  extrapolations  of  previously  established  eigenvalues  for 
fixed  R* s  and  ac*s.  The  third  is  obtained  from  a  linear  prediction  based 
on  the  first  two.  Thereafter,  the  aforementioned  procedure  is  used.  The  value 
for  a  in  the  next  trial  is  simply  a^.  It  is  given  by  the  equation 

(zrzi  -zrZiW+/zrZi  -zPz1W  +  /zrzi  -  Zr  Zi  y 

a  -  \  ?  3  3  2/  V  3  1 _ 1  ±1  V  1  2  2  V  -  (10) 

0  /zrZi  -Zrzi\  +  /zpzi  -zrz"\+/zrZi  -zxzl\ 

{  2  3  3  2)  {  }  1  1  3/  \  1  1  2  V 

Ibder  the  present  scheme,  there  is  one  free  parameter  that  one  may  fix  to 
modify  the  topology  of  the  surfaces  expressed  by  equation  (17).  That  param¬ 
eter  is  D  in  equation  (l^).  It  was  found  that  setting  D  equal  to  IT 
led  to  satisfactory  results  for  all  the  cases  considered.  It  is  sometimes 
desirable  to  select  a  different  set  of  arguments  to  solve  the  characteristic 
problem.  To  determine  points  on  a  curve  of  constant  temporal  amplification 
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such  as,  say,  the  neutral  curve  (c^  =  0),  one  may  choose  «r  and  A,  where 
A  is  defined  by  the  relation  A  =  ^r/cr*  This  procedure  was  in  fact  used  to 
generate  eigenfunction  solutions  for  neutral  disturbances. 

6.6  De termination  of  the  Eigenfunction 

Once  the  eigenvalues  are  found,  the  task  of  constructing  the  corresponding 
eigenfunction  is  relatively  simple.  The  two  solutions  denoted  by  ^  ^(y) 
and  corresponding  to  some  yfc  are  related  to  the  previous  sets 

obtained  by  the  Gram-Schmidt  orthogonalization  procedure  by  the  recursion 
formulas 


and 


~  *3,k-l  *k  S,k-2  ^  \-l 

(20a) 

m-1 

■  n  ■ 8 1 

(20b) 

j=o 

*l,k  =  *l,k-l  “  \  *3,k 

(21a) 

=  *l,k-2  “  \  “  \-l  *3,k-l 

(21b) 

m-1 

r — i 

°  *l,k-m  “  L  *k-n  *3,k-n  1 

(21c) 

n=0 

Hie  constant  is  the  normalization  factor  for  k  at  In  the 

course  of  integration,  «  (y)  and  4>,  ,  (y)  are  evaluated  at  a  single 

k-h  jy  K-m  _ 

point.  Biese  are  in  fact  the  quantities  $1(mh)  and  <|>^(mh)  defined  earlier. 
Bierefore,  it  is  convenient  to  express  equation  (21c)  in  a  slightly  different 
form.  From  equation  (20b)  it  may  be  shown  that 

T  m”1 

=  J]  m  —  1  >  n  (22) 

^3,  k-m  j=n 


IT 


so  that 


^l,k-m  \,m  ^3>k-m 


where 


I  V. 

n=0 


ra-1 

n 

d=n 


(23) 


(24) 


Note  that 

bk,m+l  =  \-m  +  \-ra-' 

The  solutions  $  (y)  and  (y)  may  be  readily  evaluated  from  the  above 

formulas  by  simply  saving  the  values  of  <j>^, 

factor  at  each  integration  step.  The  eigenfunction  is  determined  by 
combining  equations  (l3),  (20b),  and  (23)  to  yield 


V 


and  the  normalization 


where 


$>(mh)  =  (mh )  + 


m-1 

b  n  ^  -  v 

.  d=° 


<(>5(mh) 


m-1  m-1 

v  =  1 A(nh)  n 

n=0  j=n 


(26) 


(27) 


6.7  Solution  of  Velocity  Profile  of  the  Mean  Flow 

The  P&lkner-Skan  profiles  are  obtained  by  a  numerical  integration  of  the 
boundary -layer  equation 


F'"(t))  =  -  FF”  +  P(F’2  -  l) 


(28) 
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with  the  appropriate  value  of  0.  The  boundary  conditions  for  thi3  equation 
are 

F(0)  =  F*<0)  =  Lim  (F*  -  l)  =  0  (29) 

q  —  oo 

A  family  of  solutions  of  equation  (28)  can  be  generated  by  using  the  initial 
conditions  given  by  equation  (29)  and  various  trial  values  of  F*'(o).  These 
solutions  may  be  parametrically  represented  as  functions  of  r\  and  I*'(o). 
Making  use  of  the  asymptotic  nature  of  the  desired  solution  reduces  the 
original  two-point  boundary-value  problem  to  finding  the  minimum  (zero)  of 
a  function  fl  given  by 

fl[F"(0),  T|]  =  W^F’  -  l)2  +  W2(F")2  +  W^F'")2  (30) 

The  quantities  W^,  W^,  and  are  weight  factors.  Since  the  boundary 
condition  on  F*  (t^)  is  explicit,  whereas  those  on  F"(t^)  and  F,m(t^)  are 
implicit,  was  taken  as  10  and  W2  and  were  taken  as  unity.  The 
particular  form  of  the  constraint  relation  given  by  equation  (30)  was  chosen 
because  it  is  suited  for  solving  boundary-layer  problems  for  which  F*  is  a 
monotonically  increasing  function.  The  process  is  initiated  by  using  trial 
values  of  F**(o)  until  at  least  one  high  and  one  low  value  have  been  obtained. 
A  high  or  low  value  is  determined  by  integrating  equation  (28)  until  the  pro¬ 
duct  (F*  —  l)F**  becomes  positive  and  then  checking  the  sign  of  either  factor. 
If  it  is  positive,  the  value  of  F**(o)  is  high;  if  it  is  negative,  the  value 
of  F"(o)  is  low.  Quadratic  interpolation  is  employed  to  yield  G|njJ1tFM(o)] 
equal  to  zero  from  the  last  high  and  low  trials  of  ^(o).  Rapid  convergence 
of  this  process  has  been  established. 

With  F(o),  F'(o),  and  F"(o)  specified,  the  solution  of  equation  (28) 

is  treated  as  an  initial-value  problem  that  is  solved  by  a  Runge -Kutta  method 

for  third-order  differential  equations.  The  calculations  were  made  on  the 

709**  IbM  computer  with  double -precision  arithmetic.  The  profiles  obtained  are 

Ilk) 

in  agreement  with  Smith's'  ,  which  are  tabulated  to  five  decimal  places. 
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The  normalized  velocity  and  its  second  derivative,  which  appear  in  the  0-S 
equation,  are  related  to  the  function  F  as  follows: 


U(y)  =  F' (y)  (31) 

and 

tf'(y)  =  (n6)2  F"'(y)  (32) 

where  %  is  taken  to  correspond  to  the  point  at  which  the  normalized 
velocity  has  the  value  0.9990.  The  velocity  profiles  thus  obtained  are 
used  to  integrate  the  0-S  equation  across  the  boundary  layer. 

The  values  obtained  for  the  momentum  thickness  0,  shape  factor  H, 
boundary- layer  thickness  t^,  and  the  displacement  thickness  A*,  are  given 
in  table  II.  The  characteristic  length  used  in  normalizing  the  stability 
parameters  is  5*  =  A*  \/(2  —  (3)vx^/Ue . 


Table  II 

BOUNDARY -IAYER  PARAMETERS  FOR  VARIOUS  0's 


p 

00 

A*  =  J  F'dn 

0 

00 

0  =  /F'(i-F’)dn 

0 

H 

1.0 

3.1**  3 

0.6479 

0.2923 

2.216 

0.8 

3.280 

0. 6987 

0.3119 

2.240 

0.6 

3.  Mo 

0.7640 

0.3359 

2.274 

0.9 

3.533 

0.80*46 

0.3503 

2.297 

o.*4 

3.636 

O.8526 

0.3667 

2.325 

0.3 

3.752 

0.9110 

0.3857 

2.362 

0.2 

3.887 

0.9842 

0.4g82 

2.411 

0.10 

4.o48 

1.0803 

0.4355 

2.481 

0.05 

4.145 

1.1417 

0.4515 

2.529 

0.0 

4.257 

1.2168 

0.4696 

2.591 

-0.05 

4.390 

1.3124 

0.4905 

2.676 

-0.10 

4.561 

1.4427 

0.5150 

2.801 

-O.ll* 

4.744 

1.5459 

0.5386 

2.963 

-0.1988 

5.562 

2-359 

0.5854 

4.029 
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7.0  RESULTS  AND  DISCUSSION 


The  results  of  the'  present  study  are  presented  in  the  form  of  curves  of 
constant  on  the  diagrams  of  dimensionless  frequency  vs.  Reynolds  num¬ 

ber  (o>r,  R^)  in  figures  1  through  l4.  The  critical  Reynolds  number  Rg#^  ^ 
determined  by  cross -plotting  the  computed  data  is  presented  as  a  function  of 
Hartree  p's  in  table  III  and  is  shown  in  figure  15.  The  Reynolds  number  at 
maximum  spatial  amplification  rate  is  presented  as  a  function  of  0  in  fig¬ 
ure  16.  For  purposes  of  comparison,  Pretsch's  results  (obtained  on  the  basis 
of  temporal  considerations)  are  also  given.  For  the  separation  profile, 
Pretsch's  results  indicate  that  unstable  disturbances  are  always  present,  what¬ 
ever  the  Reynolds  number.  This  strongly  suggests  that  all  separation  flows  are 
essentially  turbulent.  A  lower  mode,  perhaps  corresponding  to  Pretsch's  was  in 
fact  discovered  in  the  present  analyses.  This  is  shown  in  figure  31.  The  wave¬ 
number  spectrum  for  this  mode  is  especially  small,  indicating  that  insofar  as 


Table  III 

CRITICAL  REYNOLDS  NUMBER  FOR  VARIOUS  0' s 


0 

R&*  .. 
cnt 

Present  Analysis 

Pretsch 

1.0 

12490 

12600 

0.8 

10920 

0.6 

8890 

8300 

0.5 

7680 

0.4 

6230 

5000 

0.3 

4550 

0.2 

2830 

3200 

0.1 

1380 

0.05 

865 

0.0 

520 

660 

-0.05 

318 

-0.10 

199 

126 

-0.l4 

138 

-0.1988 

67 

0 
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stability  is  concerned  the  dominant  mode  is  the  one  possessing  the  critical 
Reynolds  number  =  67  shown  on  the  same  graph.  Unstable  disturbances  cor¬ 
responding  to  the  lower  mode  were  found  to  propagate  upstream  (c  <  0)!  These 

6 

self -excited  oscillations  apparently  do  not  lead  to  transition,  at  least  until 
a  critical  Reynolds  number  is  reached,  for  laminar  separation  has  certainly 
been,  experimentally  observed.  The  maximum  frequency  for  which  disturb¬ 

ances  are  unstable,  as  well  as  the  maximum  amplification  rate  (’ai^*/R5*)max> 
is  shown  in  figure  17* 

In  order  to  carry  out  a  more  thorough  comparison  with  the  results  of 
other  investigators,  and  particularly  with  those  of  Pretsch,  curves  of  constant 
c i  and  cr  calculated  from  the  present  spatial  results  are  presented  on  the 
diagram  of  dimensionless  wave  number  vs.  Reynolds  number  (a^b*,  Rg#)>  figures 
18  through  31 •  The  neutral  curves  are  summarized  in  figure  32.  Contours  of 
the  neutral  curves  obtained  by  Pretsch  (not  shown)  are  generally  in  fair  agree¬ 
ment  with  the  present  results  for  positive  P’s.  Pretsch* s  neutral  curve  for 
the  Blasius  profile  encloses  a  much  smaller  wave  number  spectrum  for  Reynolds 
numbers  near  the  critical  as  may  be  seen  in  figure  33,  in  which  the  predicted 
maximum  wave  number  for  which  disturbances  are  unstable  (a  6*)  is  plotted 
as  a  function  of  p.  The  dontour  along  the  lower  branch  of  the  neutral  curve 
was  found  to  be  in  good  agreement  with  the  present  results.  The  maximum  tem¬ 
poral  amplification  rate,  ci]nax>  ^  a3-so  compared  in  figure  33«  Pretsch* s 
results  were  obtained  on  the  basis  of  a  temporal  criterion  and  the  present 
results  on  the  basis  of  a  spatial  criterion.  The  agreement  shown  between  the 
results  for  maximum  amplification  rates  obtained  from  the  two  analyses  is  an 
indication  that  the  use  of  the  Gaster  transformation  in  connecting  the  results 
of  the  spatial  and  temporal  analyses  is  valid,  at  least  insofar  as  unstable 
disturbances  are  concerned.  The  discrepancy  between  the  results  of  the  two 
analyses  is  perhaps  best  explained  by  the  fact  that  asymptotic  methods  are  less 
accurate  in  regions  where  viscous  effects  are  important.  In  this  region  the 
in  viscid  solution  is  less  representative  of  the  total  solution.  The  viscous 
correction  is  therefore  less  able  to  compensate  for  this  defect,  particularly 
at  low  Reynolds  number  such  as  near  the  critical  Reynolds  number  of  a  Blasius 
profile  or  profiles  having  an  adverse  pressure  gradient. 

In  the  case  of  the  negative  Hartree  P*s,  it  has  been  theoretically 
established  that  at  large  Reynolds  numbers  the  upper  branch  of  the  neutral 
curve  does  not  approach  the  R^-axis,  as  it  does  with  positive  0*s,  but 
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approaches  a  certain  asymptote  a  6*  =  constant.  This  result  was  discussed 
(l*5)  r 

by  Tollmienv  ,  who  showed  that  the  constant  is  a  function  of  P.  The 
asymptotes  for  Ot^b*  for  these  cases  provide  a  further  check  on  the  validity 
of  the  present  results,  because  the  limiting  inviscid  solution  must  emerge  as 
the  exact  solution  of  the  complete  Orr-Soramerfeld  equation  as  R  -*  «>  for 
nonzero  values  of  the  wave  number  a.  In  this  region  the  asymptotic  method 
used  by  Pretsch  should  be  valid.  The  comparison  in  figure  34  for  p  =  —  0.10 
shows  good  agreement  with  the  asymptotic  value  obtained  by  Pretsch.  It  may  be 
noted  from  figure  28  that  the  neutral  curve  for  p  =  —  0.05  obtained  from  the 
present  analysis  had  not  reached  an  asymptote  even  at  R^  =  10  .  For  this 
reason,  the  calculations  were  extended  to  a  Reynolds  number  of  200,000  for 
this  case.  Similar  calculations  were  made  for  the  cases  p  =  —  0.01  and 
P  =  —  0.025.  In  none  of  these  cases  was  an  asymptote  firmly  established, 
although  the  gradient  present  was  found  to  be  small.  The  values  for 
V  =  2  x  105  are  plotted  on  the  figure. 

Velocity  perturbations  for  the  Blasius  profile  at  Rg*  =  902  and  2080 
were  obtained  by  matching  the  minimal  value  of  the  eigenfunction  solutions  to 
the  experimentally  observed  amplitudes  of  Schubauer  and  Skramstad.  These  are 
shown  in  figures  55  and  36,  respectively.  From  the  description  given  in 
Schubauer  and  Skramstad* s  report,  one  cannot  precisely  determine  the  location 
in  the  profile  where  the  measurements  were  made.  To  locate  a  neutral  disturb¬ 
ance,  the  hot  wire  would  have  had  to  be  moved  several  inches  downstream  in 
order  to  get  a  measurable  change  in  amplitude.  As  a  result,  the  wire  would 
have  had  to  be  moved  upward  enough  to  maintain  the  same  similar  distance  from 
the  surface.  This  adjustment  is  physically  difficult  to  accomplish.  The  point 
at  which  phase  reversal  takes  place  is  more  readily  established  experimentally. 
This  may  be  done  by  simply  fairing  through  the  available  data.  The  data  of 
Schubauer  and  Skramstad  is  rescaled  to  match  the  point  of  phase  reversal  with 
the  results  of  the  theoretical  analysis  (shown  on  same  figure).  It  may  be 
observed  in  figures  35  and  36  that  the  correlation  with  theory  is  significantly 
improved  in  the  two  cases  considered.  Of  particular  note  is  the  much  closer 
agreement  obtained  for  the  location  of  the  maximum  and  minimum  within  the 
boundary  layer.  The  solutions  obtained  by  Radbill^^  have  also  been  included 
in  these  figures. 


The  agreement  between  the  results  of  the  theoretical  analyses  is  quite 
good  and  is  somewhat  reassuring,  since  in  the  present  analysis  no  investiga¬ 
tion  intp  the  effect  of  the  frequency  of  application  of  the  Gram-Schmidt 
procedure  on  the  eigenfunction  solution  had  been  made.  The  computation  time 
required  to  establish  one  set  of  eigenvalues  is  roughly  20  seconds.  A  compari¬ 
son  with  Kaplan’s  result  for  =  3000  is  shown  in  figure  37. 
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8.0  CONCLUDING  STATEMENTS 


'Solutions  for  Reynolds  numbers  up  to  Rg#  =  10^  were  successfully 
obtained  by  the  present  method  for  all  the  cases.  This  Reynolds  number  is  at 
least  50  times  as  large  as  that  which  can  be  handled  by  the  "purification” 
scheme  devised  by  Landahl  and  Kaplan.  The  method  is  capable  of  being  extended 
into  the  inviscid  region  and  in  fact  precludes  the  necessity  for  performing 
an  inviscid  analysis  in  most  practical  situations.  It  is  clear  that  the 
present  method's  capability  of  coping  with  the  numerical  difficulties 
encountered  in  performing  a  stability  analysis  of  small  disturbances  in  an 
inviscid  problem  may  be  adapted  to  deal  with  compressible  flows. 
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Figure  18.  -  Curves  of  constant  temporal  amplification  rates  (/?  =  1.0) 


Figure  21.  -  Curves  of  constant  temporal  amplification  rates  (0  -  0.5) 
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Figure  24.  -  Curves  of  constant  temporal  amplification  rates  (0  =  0.2) 


Figure  25.  -  Curves  of  constant  temporal  amplification  rates  (0  =  0.1) 
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square  longitudinal  disturbance  velocity  with  data  of  Schubauer  and  Skramstad  (R$*  ^  902) 


180*  PHASE  SHIFT  u'/Ue,  PERCENT 
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Figure  36.  -  Comparison  of  root-mean-souare  longitudinal  disturbance  velocity 
with  data  of  Schubauer  and  Skramstad  (R§*  =  2080) 
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PRESENT  METHOD 
KAPLAN  (1964) 


Figure  37.  -  Comparison  of  eigenfunction  solutions  with  the  results  of  Kaplan  (0  =  0,  R$  *  3000) 
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